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Taylor series integration is implemented in a spacecraft trajectory analysis code - the 
Spacecraft N-body Analysis Program (SNAP) - and compared with the code’s existing 
eighth-order Runge-Kutta Fehlberg time integration scheme. Nine trajectory problems, 
including near Earth, lunar. Mars and Europa missions, are analyzed. Head-to-head 
comparison at five different error tolerances shows that, on average, Taylor series is faster 
than Runge-Kutta Fehlberg by a factor of 15.8. Results further show that Taylor series has 
superior convergence properties. Taylor series integration proves that it can provide rapid, 
highly accurate solutions to spacecraft trajectory problems. 
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Nomenclature 

x component of spacecraft position relative to inertial frame centered at the central body 
y component of spacecraft position relative to inertial frame centered at the central body 
z component of spacecraft position relative to inertial frame centered at the central body 
x component of spacecraft velocity relative to inertial frame centered at the central body 
y component of spacecraft velocity relative to inertial frame centered at the central body 
z component of spacecraft velocity relative to inertial frame centered at the central body 
spacecraft mass 

(x t , X 2 , X 3 , X 4 , X 5 , X 6 , X 7 ) = spacecraft state vector 
(Xj,X 2 ,X 3 ) = spacecraft position vector 

(x 4 , X 5 , X 6 ) = (Vj , V 7 , V 3 ) = spacecraft velocity vector 
time 

d/dt = derivative with respect to time 
d/dt = derivative with respect to time 
mass flow rate 
thrust magnitude 

x component of spacecraft acceleration relative to inertial frame centered at the central body 
y component of spacecraft acceleration relative to inertial frame centered at the central body 
z component of spacecraft acceleration relative to inertial frame centered at the central body 
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= (a l ,dT,a i ) = acceleration vector 

= position vector of the f h other body relative to the central body 

= gravitational constant 
= body mass 
= step size 
= local error tolerance 
= step multiplication factor 


I. Introduction 

T he advantages of Taylor series integration in solving ordinary differential equations have been known for some 
time 1_25 . Foremost among these is the ability to maintain very high computational efficiency while achieving 
high accuracy. In fact, comparisons with other methods have shown that Taylor series integration can be faster by a 
factor of twenty or more 20 . 

The success of the method depends on recasting the governing differential system into a canonical form whereby 
system derivatives can be obtained to arbitrary order through recursion. This obviates the need to directly calculate 
derivatives, and makes it possible to obtain derivative information cheaply. The system variables can thus be 
expanded in a highly accurate series at each time level at minimal cost. 

It turns out that most differential systems can be recast into the required canonical form in a straightforward 
manner. This has led a number of authors to develop general purpose software which can recast an arbitrary system 
into canonical form and solve the resulting equations automatically Taylor series integration can thus 

be used as both a general purpose solver and also for specific applications. 

The focus here is on calculating spacecraft trajectories. Previous work using Taylor series to calculate 
trajectories includes that in Refs. 21 and 22. Unlike Refs. 21 and 22, which used an automated Taylor series 
package 16 22 to propagate trajectories in Earth orbit, the present work focuses on using Taylor series integration in 
an existing trajectory analysis code, SNAP (Spacecraft N-body Analysis Program) 26 . Developed at NASA’s Glenn 
Research Center, SNAP is a high fidelity trajectory propagation program that can propagate the trajectory of a 
spacecraft about virtually any body in the solar system. The equations of motion include the effects of central body 
gravitation with N x N harmonics, other body gravitation with N x N harmonics, solar radiation pressure, 
atmospheric drag (for Earth orbits) and spacecraft thrusting (including shadowing). The equations are solved using 
an eighth-order Runge-Kutta Fehlberg (RKF) 27 single step method with variable step size control. 

The purpose of this paper is to demonstrate the use of Taylor series (TS) integration in a high fidelity trajectory 
analysis code, SNAP, and to provide a detailed comparison of TS performance to eighth-order RKF ~ 7 . Section II 
presents the equations of motion, Section III describes the TS formulation and Section IV discusses the numerical 
implementation. Section V compares TS and RKF on a representative set of spacecraft trajectory problems, 
including near Earth, lunar, Mars and Europa missions. It is shown that TS is faster than RKF by an average factor 
of 15.8, while simultaneously improving accuracy. 


II. Equations of Motion 

Let X = (x l ,x 2 ,x 3 ,x 4 ,x 5 ,x 6 ,x 7 ) denote the spacecraft state vector, where (x, , X 2 , x 3 ) = x is the 
spacecraft position in Cartesian coordinates relative to an inertial frame centered at the central body, 
(x 4 , X 5 , X 6 ) = V is the spacecraft velocity relative to an inertial frame centered at the central body, and X 7 is the 
spacecraft mass. The equations of motion are 


Xj = x 4 

x' = x 5 (1) 

X3 = v 6 
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where fl ( is the acceleration in the / th coordinate direction and til is the mass flow rate. The acceleration is a 

function of the forces acting on the spacecraft. Forces included in SNAP are central body, other body, thrust, 
atmospheric drag (for low Earth orbits), solar radiation pressure, oblateness effects of Earth, and oblateness effects 
of other bodies, so that 


a = a cb + a ob + a,h + a d + a srp + a obe + a obo 


( 2 ) 


This paper considers only the first four acceleration terms, which are given below. 


a chJ = - G M cb 


x. 


(x 2 + x 2 + x 2 ) 2 


T : = ~ GM cbJzk 

X 


( 3 ) 


where G is the gravitational constant, M cb is the central body mass, i = 1 ,2,3 denotes the coordinate direction, and 
the spacecraft mass is much less than the central body mass. 


\b,t=^ J - GM j 


l 


X i~ X U , X iJ 

3 + _ 3 
X - X , X. 


\ 


where j denotes the / h body and lx 


: = ( x lj +x h +x lj ) 2 - 


T v i 1 

a th,i =r R — 

v x. 


where T is the constant thrust magnitude and the thrust direction is parallel to the velocity vector. 

a d\ = C 1 |( X 4 + C 2 X 2 y + ( X 5 - C 2 X 1 y + X 6 J ( X 4 + C 2 X 2 V X 7 
a d,2 = C 1 |( X 4 + C 2 X 2 f + ( X 5 - C 2 X 1 y + X 6 J ( X 5 “ C 2 X 1 V X 7 

3 = C 1 |( X 4 + C 2 X 2 y + ( X 5 - C 2 X 1 y +X l\ X J X 7 


( 4 ) 


( 5 ) 


( 6 ) 


( 7 ) 


( 8 ) 


where C, = - — C 7) • ( spacecraft area ) • ( atmospheric density ) , c 2 = rotation rate of Earth and C w is the 


drag coefficient. 
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III. Taylor Series Formulation 


Let the state vector X have initial condition X 0 . Within the radius of convergence, the system variables 
X n ( t ) can be expanded in a Taylor series, 

" v-WfV ) 

*„0> = X " k\ ^~ to ^ k ’ n = l ’-’ 7 (9) 

(k) 

where the derivatives X n are obtained by successively differentiating the right hand side of Eqs. (1). This can be 

efficiently accomplished using recurrence relations, as follows. Consider only the central body acceleration term, so 
that 


Introduce X 8 and x g . 


Eqs. (10) become 


(10) 


(ID 


( 12 ) 


( 13 ) 


NASA/TM— 2008-2 1 5439 


4 



x ' 6 =-GM cb ^ 

x 9 

x'i = 0 


and two auxiliary equations are added to the system: 


Xg = 2x x x 4 + 2x 2 x 5 + 2x 3 x 6 


f 3 XgXg 


x 9 = 


2 x c 


(14) 

(15) 


The right hand side of the new system, Eqs. (13) - (15), can now be differentiated using recurrence relations for 
products and quotients. 

For a function w{t) = f ( t ) g{t ) , the Leibnitz rule for differentiating products gives 16 


W{k) = ^ j F(j)G{k-j) 


j = o 


(16) 


W W (( 0 ) 


where W (k) := 


k\ 


FUY-F1M and 


J 




are reduced derivatives. 


fit) 


For a function w(t) = J - ’ , the recurrence relation for quotients gives 1 1 

gf) 


W(k) = - 

g 


F(k)-^G(jW(k~j) 

j-i 


(17) 


where W, F and G are reduced derivatives as above. 
The recurrence relations are derived as follows. Let 


Eqs. (13) — (15) become 


x', n = 

(18) 

X 

ii 

s" 

(19) 

u 2 = x 5 

(20) 

m 3 = x 6 

(21) 

u a --G M cb — 
x 9 

(22) 

u 5 = ~G M cb — 
x 9 

(23) 
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( 24 ) 


u 6 = ~G M cb — 

Xg 


u 1 = 0 

(25) 

u 8 = 2xjX 4 + 2x,x 5 + 2x 3 x 6 

(26) 

3 XgVl g 

2 x 8 

(27) 

X, X 7 X, 

Introduce auxiliary functions W 4 = — , W 5 = — , W 6 = — , W S1 = X jX 4 , 

Xg Xg Xg 

XqUo 

andWg = 

x 8 

Using Eq. (18) one obtains 

w 8 , 2 = x 2 x 5 , w 8>3 = x 3 x 6 , 

«r } =x^ 

(28) 

from which 


w (iM) x w 

" = " , k > 1 

(A'-l)! (A: - 1)! 

(29) 

=► u„(k-\)-kxsk) =» r.(t)-EdUd) , 

k 

k > 1 (30) 

The recurrence relations are then, for all k > 1 , 


u l (k)-x i (k)- u * {k ~ l) 

k 

(31) 

k 

(32) 

U 3 (k)~X,(k)~ Ui(k ~ l) 

k 

(33) 

U 4 (k) = - G M cb W 4 (k) 

(34) 

U 5 (k) = - G M cb W 5 (k) 

(35) 

U 6 (k) = - G M cb W 6 (k) 

(36) 
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U 1 (k) = 0 


(37) 


where 


U t (k) = 2 W sl (k) + 2 W St2 (k) + 2 W 8 3 (k) 
U 9 (k) = jw 9 (k) 


(38) 

(39) 


W A (k) = — 

Xq 


X,{k)-^X 9 (j)W A {k- j) 

M 


Xa 


u x {k- 1) At/ 9 (i-l) 


-2 


W-y) 


ft. 


m~^x,(j)x A (k-j) = 


C/ 4 (A:-1) 


x 4 + 


M J 

v^iO'- 1) u A {k-j-\) 


7=o 


-i 7 


k~j 


(40) 


(41) 


etc., and a similar expression can be derived for IF 9 . 
The Taylor series coefficients are then 




7 r >*„(*) - 

k\ 


U„ (k - 1) 
k 


l< k < K 


and the local series solution is 


*„(0 = 


A- „(*■) 


K Vo ) 

k\ 


it-t 0 ) k + ^ 


(42) 


(43) 


where U n (k) is defined by Eqs. (31)- (39), U n (0) is defined by Eqs. (19) - (27), K is the number of terms in the 
series and T n K is the truncation error. 

The other acceleration terms can be handled similarly. Only other body acceleration, Eq. (4), requires special 
consideration, due to the need for the motion of other bodies. This can generally be obtained from ephemeris files. 
However, integration by Taylor series requires derivatives not available from ephemeris files. It is thus necessary to 
integrate the other body motion as part of the governing differential system. This leads to a substantially larger 
system of equations, but fortunately can still be integrated efficiently. 

Once the recurrence relations are derived for all acceleration terms and the state vector specified, Eqs. (42) - 

(43) are used to expand the system variables in a series from t 0 to t t , where the step size h x := — t 0 is 

determined to meet the local error tolerance. From t x , the variables are expanded in a new series to t 2 , and so 
forth. Thus, by a process of “analytic continuation,” one obtains a set of overlapping series solutions that cover the 
integration domain. 


IV. Numerical Implementation 

Taylor series integration was implemented in SNAP by making some minor modifications to existing source 
code and adding three additional subroutines - a driver routine which automatically introduces auxiliary variables, 
sets up initial conditions and integrates; a routine which calculates system reduced derivatives using recurrence 
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relations for the following quotients and products: — — , X m X n , m n , X m x' m , — — 9 m n and X m X n ; and a 

x x xx, 

n n n l 

routine which determines the step size and sums the series. The number of series terms is variable up to a maximum 
of 30, but remains constant throughout the integration. Positive and negative terms are summed separately to avoid 
cancellation of significant digits. 

The step size can be determined from the standard formula 28 


^nex, = r l h 


\ M 


(44) 


where h denotes the current step, T the local error tolerance, e max the estimate of maximum truncation error, M 

the order of the maximum truncation error estimate and Tj < 1 the step multiplication factor. Eq. (44) is more or 

less restrictive depending on I] and the truncation error estimate e max . Generally, e max should not be calculated 

from the next series term, due to the extra computation required and the fact that it is not a reliable error estimate 29 . 
A conservative approach which takes advantage of the series terms already computed leads to 


- Max, \_\X„(K-V,\h K -' + ly.ml/i"] 


(45) 


where the expression in brackets is derived by subtracting the Taylor series solution of degree K - 2 from the 
solution of degree K and taking absolute values of individual terms. Eq. (45) can be viewed as a truncation error 
estimate for the series of degree K—2 which is then applied to the more accurate series of degree K. 

An alternative to Eqs. (44) - (45) is to simply require h to be small enough that the system variables directly 
satisfy the absolute error tolerance requirement 

\X„{K-l)\h K - x + \X n {K)\h K <x (46) 


for all n. Eq. (46) can be solved by fixed point iteration, 


K i = ex P 


( i 


In, 


K- 1 X n (K - 1 ) + hj X n (K) 


(47) 


The smallest h is chosen over all n, and multiplied by the step multiplication factor t) . This approach offers the 
advantage of directly calculating the step size without the need for a previous step, and guarantees that the error 
tolerance is met. Eq. (44), on the other hand, requires a previous step and will require a repeat step whenever 

e > T . 

max 

The step selection methods above performed very similarly in the current study. Both methods provided stable, 
accurate solutions and used approximately the same number of time steps in head-to-head calculations. 


V. Results 

We compare RKF and TS performance on the trajectory problems in Table 1 . All calculations were run on a Dell 
PowerEdge 2600 with two 3.066 GHz processors and four GB of RAM. Source code was compiled using the 
Absoft Fortran 90 compiler without optimization. SNAP was run with all intermediate print and stop options turned 
off. All TS calculations used a series with 20 terms and a variable step size determined by Eq. (47) with a step 
multiplication factor that ranged from 0.75 to 0.9. 
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Table 1 


Problem 

Title 

Description 

Central 

Body 

Other Bodies 

1 

Satellite in low 
Earth orbit 
(LEO) 

A 10,000 kg satellite orbits Earth for 10 days 
at an inclination of 28.45 degrees. 

Earth 

Moon 

2 

Satellite in LEO 
with drag 

A 10,000 kg satellite orbits Earth for 10 days 
with constant drag at an inclination of 28.45 
degrees. 

Earth 

Moon 

3 

Spacecraft 
spiraling out of 
Earth’s gravity 

well 

A 10,000 kg spacecraft spirals out of Earth’s 
gravity well in a low thrust trajectory. 
Calculation stops when the semi-major axis 
of trajectory equals 40,000 km. 

Earth 

Sun, Moon 

4 

Spacecraft from 
near Earth to 
lunar orbit 

A 3580 kg spacecraft 400 km above Earth has 
been propelled with sufficient energy to reach 
the Moon. Spacecraft coasts to Moon, 
performs insertion bum, propagates to 
apolune, and performs final bum to achieve 
500 km by 10,000 km polar lunar orbit with 
an argument of perilune equal to 90 degrees. 
See Fig. 1. 

Moon 

Earth, Sun 

5 

Spacecraft in 
lunar orbit 

Spacecraft with 2848.56 kg mass coasts for 
10 days in 500 km by 10,000 km polar lunar 
orbit with an argument of perilune equal to 90 
degrees. See Fig. 1. 

Moon 

Earth, 

Sun 

6 

Spacecraft 
thrusting from 
near Earth to 
Mars coast 

A 585 kg spacecraft near Earth thrusts for 
38.45 days to achieve sufficient energy to 
coast to Mars. See Fig. 2. 

Sun 

Earth, Moon, Venus, Mars, 
Jupiter barycenter, 
Saturn barycenter 

7 

Spacecraft coast 
to Mars flyby 

A 555.66 kg spacecraft coasts to Mars flyby 
for 161.55 days. See Fig. 2. 

Sun 

Earth, Moon, Venus, Mars, 
Jupiter barycenter, 
Saturn barycenter 

8 

Spacecraft 
thrusting 
tangentially out 
of Europa orbit 

A 10,000 kg spacecraft in Europa orbit 
thmsts tangentially to spiral out until the 
semi-major axis equals 10,000 km. 

Europa 

Jupiter, Sun, Ganymede, Io 
Callisto 

9 

Spacecraft coast 
near Europa 

A 9800.49 kg spacecraft coasts for one day 
after spiraling out of Europa orbit. 

Europa 

Jupiter, Sun, Ganymede, Io 
Callisto 




Figure 1. Earth to Moon Trajectory Figure 2. Earth to Mars Flyby 
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Each trajectory was integrated at five error tolerances from 10 10 to 10 14 . Tables 2-10 summarize the results. 
RKF results are shown on top and TS on bottom. Spacecraft positions are in kilometers and CPU ratio is RKF/TS. 
RKF and TS velocities agreed equally well as spacecraft position, but are omitted here for brevity. 


Table 2 Results for Problem 1 


Spacecraft Position at End of Calculation 

CPU 

CPU 

T x 

y 

z 

(sec) 

ratio 

l.E-10 -4505.4956174253 

4460.8015923075 

2416.2243468200 

6.91 


l.E-11 -4505.4976241104 

4460.8000212437 

2416.2234952381 

9.10 


l.E-12 -4505.4977855238 

4460.7998948715 

2416.2234267390 

12.23 


l.E-13 -4505.4977981582 

4460.7998849798 

2416.2234213773 

16.05 


l.E-14 -4505.4977991371 

4460.7998842135 

2416.2234209620 

21.17 


l.E-10 -4505.4893402433 

4460.8066582279 

2416.2273550484 

0.20 

34.6 

l.E-11 -4505.4893401709 

4460.8066582847 

2416.2273550792 

0.22 

41.4 

l.E-12 -4505.4893402598 

4460.8066582151 

2416.2273550415 

0.25 

48.9 

l.E-13 -4505.4893403875 

4460.8066581153 

2416.2273549873 

0.28 

57.3 

l.E-14 -4505.4893402869 

4460.8066581940 

2416.2273550300 

0.31 

68.3 


Table 3 Results for Problem 2 

Spacecraft Position at End of Calculation 

CPU 

CPU 

T x 

y 

z 

(sec) 

ratio 

l.E-10 -6286.5348347365 

2234.8075149479 

1209.8215168263 

7.04 


l.E-11 -6286.5358380683 

2234.8053241518 

1209.8203296255 

9.16 


l.E-12 -6286.5359189241 

2234.8051476038 

1209.8202339534 

12.18 


l.E-13 -6286.5359251975 

2234.8051339059 

1209.8202265305 

16.04 


l.E-14 -6286.5359257481 

2234.8051327038 

1209.8202258791 

21.28 


l.E-10 -6286.5317837499 

2234.8146866042 

1209.8256504332 

0.26 

27.1 

l.E-11 -6286.5317836592 

2234.8146868021 

1209.8256505404 

0.29 

31.6 

l.E-12 -6286.5317837023 

2234.8146867081 

1209.8256504895 

0.33 

36.9 

l.E-13 -6286.5317837303 

2234.8146866471 

1209.8256504564 

0.37 

43.4 

l.E-14 -6286.5317836427 

2234.8146868382 

1209.8256505600 

0.42 

50.7 


Table 4 Results for Problem 3 



Spacecraft Position at End of Calculation 

CPU 

CPU 

T 

X 

y 

z 

(sec) 

ratio 

l.E-10 

21783.589218926 

30426.566335107 

14117.151742798 

44.01 


l.E-11 

21783.168251826 

30426.814392403 

14117.266771932 

58.17 


l.E-12 

21783.134612346 

30426.834214343 

14117.275963762 

77.72 


l.E-13 

21783.131993250 

30426.835757635 

14117.276679417 

102.06 


l.E-14 

21783.131799858 

30426.835871591 

14117.276732261 

134.87 


l.E-10 21783.126013948 

30426.839470527 

14117.277961164 

2.29 

19.2 

l.E-11 

21783.126025067 

30426.839463974 

14117.277958125 

2.57 

22.6 

l.E-12 

21783.126012532 

30426.839471362 

14117.277961551 

2.88 

27.0 

l.E-13 

21783.126009565 

30426.839473099 

14117.277962357 

3.27 

31.2 

l.E-14 21783.126006877 

30426.839474690 

14117.277963095 

3.69 

36.6 
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Table 5 Results for Problem 4 


Spacecraft Position at End of Calculation 


T x 

l.E-10 -178.63870529793 
l.E-11 -178.63855098380 
l.E-12 -178.63855046608 
l.E-13 -178.63855601675 
l.E-14 -178.63854695266 


y 

4876.3161604664 

4876.3170369698 

4876.3170506570 

4876.3170085371 

4876.3170607540 


z 

-10650.392195425 

-10650.391796702 

-10650.391790443 

-10650.391809627 

-10650.391785873 


CPU 

(sec) 

0.44 

0.51 

0.59 

0.71 

0.87 


CPU 

ratio 


l.E-10 -178.71939471463 
l.E-11 -178.71939481530 
l.E-12 -178.71939478451 
l.E-13 -178.71939478450 
l.E-14 -178.71939478714 


4876.1033547743 

4876.1033542962 

4876.1033544409 

4876.1033544292 

4876.1033544962 


-10650.488272923 0.19 

-10650.488273141 0.20 

-10650.488273078 0.20 

-10650.488273079 0.20 
-10650.488273051 0.20 


2.31 

2.55 
2.95 

3.55 
4.35 


Table 6 Results for Problem 5 


Spacecraft Position at End of Calculation CPU CPU 


T x 

l.E-10 -215.32731201862 
l.E-11 -215.32794605537 
l.E-12 -215.32799474610 
l.E-13 -215.32799868206 
l.E-14 -215.32799901817 

l.E-10 -215.32849813488 
l.E-11 -215.32849812774 
l.E-12 -215.32849812980 
l.E-13 -215.32849813332 
l.E-14 -215.32849813453 


y 

-1650.3697144277 

1635 

-1650.3715299629 

1635 

-1650.37166939046 

1635 

-1650.3716806616 

1635 

-1650.3716816241 

1635 

-1650.3731189359 

1635 

-1650.3731189154 

1635 

-1650.3731189214 

1635 

-1650.3731189310 

1635 

-1650.3731189349 

1635 


z 

(sec) 

ratio 

1626824928 

1.99 


1612149497 

2.64 


1611022567 

3.49 


1610931478 

4.62 


1610923701 

6.15 


1600219626 

0.32 

6.22 

1600219790 

0.35 

7.54 

1600219744 

0.40 

8.73 

1600219653 

0.44 

10.5 

1600219635 

0.49 

12.6 


Table 7 Results for Problem 6 


Spacecraft Position at End of Calculation 


T x 

l.E-10 21572817.105605 
l.E-11 21572817.105556 
l.E-12 21572817.105508 
l.E-13 21572817.105509 
l.E-14 21572817.105496 


y 

-139377155.25969 

-139377155.25962 

-139377155.25961 

-139377155.25960 

-139377155.25961 


z 

-62943596.519662 

-62943596.519125 

-62943596.519004 

-62943596.518999 

-62943596.518987 


CPU 

(sec) 

0.14 

0.15 

0.21 

0.25 

0.34 


CPU 

ratio 


l.E-10 21572816.904607 
l.E-11 21572816.904606 
l.E-12 21572816.904607 
l.E-13 21572816.904607 
l.E-14 21572816.904607 


-139377156.15407 

-139377156.15407 

-139377156.15407 

-139377156.15407 

-139377156.15407 


-62943596.9802248 0.11 1.27 

-62943596.9802247 0.12 1.25 

62943596.9802248 0.14 1.50 

-62943596.9802248 0.16 1.56 

-62943596.9802247 0.17 2.00 
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Table 8 Results for Problem 7 



Spacecraft Position at End of Calculation 

CPU 

CPU 

T 

X 

y 

z 

(sec) 

ratio 

l.E-10 

174298859.53574 

117412103.45296 

49149031.016983 

0.14 


l.E-11 

174298859.53061 

117412103.45425 

49149031.017954 

0.15 


l.E-12 

174298859.53010 

117412103.45436 

49149031.018040 

0.17 


l.E-13 

174298859.53005 

117412103.45437 

49149031.018048 

0.18 


l.E-14 

174298859.53003 

117412103.45437 

49149031.018052 

0.23 


l.E-10 

174298891.93452 

117412111.58059 

49149027.497160 

0.11 

1.27 

l.E-11 

174298891.93452 

117412111.58059 

49149027.497159 

0.13 

1.15 

l.E-12 

174298891.93452 

117412111.58059 

49149027.497160 

0.15 

1.13 

l.E-13 

174298891.93453 

117412111.58059 

49149027.497159 

0.17 

1.06 

l.E-14 

174298891.93453 

117412111.58059 

49149027.497159 

0.17 

1.35 



Table 9 Results for Problem 8 




Spacecraft Position at End of Calculation 

CPU 

CPU 

T 

X 

y 

z 

(sec) 

ratio 

l.E-10 

7528.3536043710 

1375.7425889427 

-383.13557482604 

0.79 


l.E-11 

7528.3536031572 

1375.7425956652 

-383.13557244841 

1.00 


l.E-12 

7528.3536030460 

1375.7425962755 

-383.13557223246 

1.33 


l.E-13 

7528.3536030359 

1375.7425963294 

-383.13557221337 

1.74 


l.E-14 

7528.3536030345 

1375.7425963361 

-383.13557221099 

2.25 


l.E-10 

7528.7471022132 

1375.7729665868 

-383.16004903160 

0.14 

5.64 

l.E-11 

7528.7471022129 

1375.7729665694 

-383.16004903738 

0.15 

6.67 

l.E-12 

7528.7471022110 

1375.7729665782 

-383.16004903433 

0.16 

8.31 

l.E-13 

7528.7471022134 

1375.7729665759 

-383.16004903524 

0.18 

9.67 

l.E-14 

7528.7471022136 

1375.7729665815 

-383.16004903345 

0.21 

10.7 


Table 10 Results for Problem 9 

Spacecraft Position at End of Calculation 

CPU 

CPU 

T x 

y 

z 

(sec) 

ratio 

l.E-10 -6179.7835988316 

19484.717456953 

11599.814016718 

0.12 


l.E-11 -6179.7835986778 

19484.717455196 

11599.814015865 

0.14 


l.E-12 -6179.7835986665 

19484.717454928 

11599.814015733 

0.17 


l.E-13 -6179.7835986627 

19484.717454902 

11599.814015720 

0.22 


l.E-14 -6179.7835986621 

19484.717454906 

11599.814015723 

0.27 


l.E-10 -6179.2994738523 

19483.004328592 

11598.975228083 

0.04 

3.00 

l.E-11 -6179.2994738523 

19483.004328592 

11598.975228083 

0.05 

2.80 

l.E-12 -6179.2994738525 

19483.004328591 

11598.975228082 

0.05 

3.40 

l.E-13 -6179.2994738524 

19483.004328591 

11598.975228083 

0.04 

5.50 

l.E-14 -6179.2994738522 

19483.004328591 

11598.975228083 

0.05 

5.40 


Table 11 summarizes the CPU ratios. TS is faster than RKF by more than an order of magnitude in 18 of 45 
cases. The average speedup is 15.8. For the interplanetary trajectory problems (4-9) the average speedup is 4.53. 
The gain here is smaller due to the additional equations that TS must solve to account for other body motion. As 
noted previously, TS must integrate other body motion as part of the differential system, whereas RKF obtains other 
body motion from ephemeris files. This difference in the integration methods explains the small differences in 
spacecraft positions observed in Tables 2-10. 
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Table 11 RKF/TS CPU ratios 

Problem 


T 

1 

2 

3 

4 

5 

6 

7 

8 

9 

L 

l.E-10 

34.6 

27.1 

19.2 

2.31 

6.22 

1.27 

1.27 

5.64 

3.00 

l.E-11 

41.4 

31.6 

22.6 

2.55 

7.54 

1.25 

1.15 

6.67 

2.80 

l.E-12 

48.9 

36.9 

27.0 

2.95 

8.73 

1.50 

1.13 

8.31 

3.40 

l.E-13 

57.3 

43.4 

31.2 

3.55 

10.5 

1.56 

1.06 

9.67 

5.50 

l.E-14 

68.3 

50.7 

36.6 

4.35 

12.6 

2.00 

1.35 

10.7 

5.40 


Another important property to consider is convergence. Tables 12-20 present the number of converged digits 

obtained for each spacecraft coordinate at each error tolerance, where the T = 10 14 case was used as the fully 
converged solution. RKF results are on top and TS on bottom. TS has more converged digits than RKF in 103 out 
of 108 cases, while RKF has more converged digits in one case. On average, TS has 2.63 more converged digits per 
case. The results also indicate that TS solutions are nearly fully converged at all error tolerances, suggesting that the 
step selection method may be too conservative. Finally, it should be noted that convergence itself does not 
necessarily imply accuracy. However, it does indicate that a necessary condition for accuracy is satisfied. 


Table 12 Number Table 13 Number Table 14 Number Table 15 Number 

of Converged Digits of Converged Digits of Converged Digits of Converged Digits 

for Problem 1 for Problem 2 for Problem 3 for Problem 4 


T 

X 

y 

z 

T 

X 

y 

z 

T 

X 

y 

z 

T 

X 

y 

z 

l.E-10 

6 

6 

6 

l.E-10 

6 

6 

6 

l.E-10 

4 

5 

5 

l.E-10 

6 

6 

8 

l.E-11 

7 

7 

7 

l.E-11 

7 

7 

7 

l.E-11 

5 

6 

6 

l.E-11 

6 

7 

9 

l.E-12 

8 

8 

8 

l.E-12 

8 

8 

9 

l.E-12 

7 

6 

7 

l.E-12 

6 

8 

10 

l.E-13 

9 

9 

10 

l.E-13 

9 

9 

9 

l.E-13 

8 

8 

9 

l.E-13 

6 

7 

9 

l.E-10 

10 

11 

11 

l.E-10 

10 

10 

9 

l.E-10 

10 

10 

10 

l.E-10 

9 

9 

11 

l.E-11 

10 

10 

10 

l.E-11 

10 

11 

10 

l.E-11 

9 

9 

10 

l.E-11 

10 

10 

12 

l.E-12 

11 

11 

11 

l.E-12 

10 

10 

9 

l.E-12 

10 

10 

10 

l.E-12 

10 

10 

12 

l.E-13 

10 

10 

11 

l.E-13 

10 

10 

9 

l.E-13 

10 

10 

10 

l.E-13 

10 

10 

12 


Table 16 Number Table 17 Number Table 18 Number Table 19 Number 

of Converged Digits of Converged Digits of Converged Digits of Converged Digits 

for Problem 5 for Problem 6 for Problem 7 f° r Problem 8 


T 

X 

y 

z 

T 

X 

y 

z 

T 

X 

y 

z 

T 

X 

y 

z 

l.E-10 

5 

6 

6 

l.E-10 

10 

12 

10 

l.E-10 

10 

li 

10 

l.E-10 

9 

8 

8 

l.E-11 

6 

7 

7 

l.E-11 

10 

13 

11 

l.E-11 

11 

12 

11 

l.E-11 

10 

10 

9 

l.E-12 

7 

8 

8 

l.E-12 

10 

14 

12 

l.E-12 

12 

13 

11 

l.E-12 

11 

11 

10 

l.E-13 

9 

9 

9 

l.E-13 

12 

13 

12 

l.E-13 

12 

14 

13 

l.E-13 

11 

11 

11 

l.E-10 

12 

11 

13 

l.E-10 

14 

14 

14 

l.E-10 

13 

14 

13 

l.E-10 

12 

11 

11 

l.E-11 

11 

11 

11 

l.E-11 

13 

14 

14 

l.E-11 

13 

14 

14 

l.E-11 

12 

11 

10 

l.E-12 

11 

11 

11 

l.E-12 

14 

14 

14 

l.E-12 

13 

14 

13 

l.E-12 

12 

12 

11 

l.E-13 

11 

12 

11 

l.E-13 

14 

14 

14 

l.E-13 

14 

14 

14 

l.E-13 

12 

12 

10 
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VI. Conclusion 

Taylor series integration was implemented in a high fidelity trajectory analysis 
code (SNAP) and compared with 8 th order Runge-Kutta Fehlberg on a representative 
set of trajectory problems. On average, TS was more than an order of magnitude 
faster than RKF. TS also showed superior convergence properties, having more 
converged digits than RKF in 103 out of 108 cases. Taylor series integration thus 
proved that it can provide rapid, highly accurate solutions to spacecraft trajectory 
problems. This is consistent with other reports which have found Taylor series 
integration to be superior to conventional methods in both speed and accuracy 1U6 ' 20 . 


Table 20 Number 
of Converged Digits 


for Problem 

9 


T 

X 

y 

z 

l.E-10 

10 

9 

10 

l.E-11 

11 

11 

11 

l.E-12 

11 

12 

12 

l.E-13 

12 

12 

13 

l.E-10 

13 

13 

14 

l.E-11 

13 

13 

14 

l.E-12 

12 

14 

13 

l.E-13 

13 

14 

14 
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